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We present a numerical scheme for calculating the first quantum corrections to the properties of 
static solitons. The technique is applicable to solitons of arbitrary shape, and may be used in 3+1 
dimensions for multiskyrmions or other complicated solitons. We report on a test computation in 
1+1 dimensions, where we accurately reproduce the analytical result with minimal numerical effort. 
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I. INTRODUCTION 

: 

o\ . 

Qv^ , The quantum fate of static solitons has long been studied, and beautiful results have emerged in exactly solvable two 
7-H ' dimensional field theories. In higher dimensions, very little is known except in certain very special supersymmetric 
theories, or in situations of spherical symmetry. In many theories, the quantum properties of the soliton-bearing 
sector are interesting, but beyond the reach of current calculations. One example of such a theory is the chiral model 
of nuclear physics, the solitonic sector of which describes nucleons and nuclei Jl]. It is clear that quantum effects must 
be included if the theory is to match real nuclei, but hitherto only a very limited 'collective coordinate' quantization 
has been possible. There are many other 3+1 dimensional field theories with topological solitons, for example those 
bearing magnetic monopoles and vortex strings. It would be very useful to develop techniques to study quantum 
corrections to these solitons in both supersymmetric and non-supersymmetric contexts. 

The quantum corrections to a soliton's classical properties may be expressed as an expansion in a dimensionless 
parameter, some power of the coupling constant, multiplied by Planck's constant h. At weak coupling, the corrections 
are small and may be calculated perturbatively in this small parameter. In this paper we outline a straightforward 
method for numerically computing the first quantum corrections, including all fluctuations about a static soliton, to 
H ' nrs t order in Ti. In particular, we present a method for calculating the quantum correction to the soliton mass, although 
the technique is easily generalized to other quantum corrections. We apply the method to a simple two dimensional 
theory where the exact expression is known, and show that it is reproduced numerically with little computational 
\ cost. We also discuss the extension of the method to analogous calculations in four dimensions. 

There exist a set of standard techniques for calculating the quantum corrections to the mass, moments, and other 
Q_i! properties of solitons 0. The techniques usually used require knowledge of the phase shifts associated with every 
possible meson scattering from the soliton. However in many cases of interest, for example multi-Skyrmions 0, the 
soliton field is not spherically symmetric and the necessary information would be extremely difficult to extract. 

In this paper we present a different approach, based upon quantizing the soliton in a finite box, which provides 
an infrared regulator and reduces the number of modes to be quantised to a discrete set. Our method centers upon 
a formula for the first quantum mass corrections due to Cahill et al. Q, which automatically removes the worst 
divergences. The quantum correction to the soliton mass is given as a trace over a complete set of modes. We 
compute this trace by summing first over the lowest normal modes of the soliton, and then over a plane wave basis up 
to a finite cutoff. The contribution of modes beyond this cutoff is included analytically using a derivative expansion. 
As we shall show, one per cent accuracy in the mass correction is achievable with very modest computational resources. 
We point toward the application of this technique to more interesting problems, where it appears feasible to calculate 
the quantum corrections for multisolitons in the Skyrme model 01, and other 3+1 dimensional theories. 



II. FLUCTUATIONS ABOUT CLASSICAL SOLITONS 



We would like to quantize small fluctuations about some stable, static classical soliton. The field's dynamics are 
determined by a Lagrangian, C(<j), d(f), which gives classical equations of motion 
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(i) 



with the soliton 4>g t (x) as a time-independent solution. The soliton is a smooth, localized 'lump' in space, surrounded 
by a vacuum of constant field. The field <fi a (x,t) may be represented as the sum of the (c-number) classical soliton, 
4> s t(x), plus a quantum correction e(x 7 t) which is an operator obeying canonical commutation relations. We work to 
first order in h, which amounts to keeping terms up to second order in e(x,t) in the the Hamiltonian. 
In this approximation, the quantum mechanical Hamiltonian is 



H T = E cl + : J d d x n 2 (n a ,e a ,d l e a ,x) : , 



(2) 



where E c \ = J d d xV(4> s t,d(f> s t) is the classical soliton mass, and TL2 quadratic in e and its conjugate momentum it. 
Normal ordering is with respect to the trivial vacuum normal modes for e, the free mesons, so the Hamiltonian is zero 
in the soliton-free vacuum. 

We may express e a (x, t) in terms of normal modes about the classical soliton. It obeys a linear equation 



H ab e 



(3) 



with H 2 b a positive semi-definite differential operator depending on the soliton solution. The eigenvalues of H 2 are 



the frequencies squared ui 2 of the normal modes about the soliton. We then write 

-iu n t 



a n e^(x) + c.c. 



(4) 



where e®(x) is the nth classical normal mode, and a n is the associated annihilation operator. In terms of these normal 
mode operators, the Hamiltonian is 



E r 



E 



w„ala 



(5) 



where normal ordering is with respect to the trivial soliton-free vacuum. The ground state of this Hamiltonian in the 
presence of a soliton can be computed by normal ordering the second term with respect to the soliton normal modes. 
To do so, write the operators a n , in terms of a a ^(k), a a (k), the creation and annihilation operators for fluctuations 
about the soliton-free vacuum. 



l n 9 E 



3(-fc)a't(fc) ( \l-r + \ ~r ) + 3(*)a°(*) (J— -J — 



(6) 



where lj^ = \/k 2 + m 2 are the frequencies of elementary meson excitations. Applying this identity, one finds 

Aa n -.+W^)el{-k) {Un ~ Wk? 



: ijj n a' n an : +-{n 



(H - H Q ) 



Ho 



(7) 
(8) 



where H is the single particle Hamiltonian for fluctuations about the soliton, i.e., the positive square root of H 2 , and 
H Q is the single particle Hamiltonian for fluctuations about the soliton-free vacuum. 
So, the total Hamiltonian may be written 



H T = E c \ + ^ Indian + 5m, 



(9) 



where 5 m , the quantum correction to the soliton mass, is 



5m = — Tr 
4 



(H - H ) 



- 1 



Hn 



(10) 



a result obtained by Cahill, Comtet and Glauber M. The trace is taken over any complete set of states. In one 
dimension, (^) is finite and no further renormalization is needed. In three space dimensions it diverges logarithmically, 
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and so requires regularization, which may be provided by inserting a factor e _ . The logarithmic divergence which 
results may then be subtracted with a local counterterm Jl2|]. Nevertheless, the formula (|l0|) is still very useful in 
three dimensions, as it has removed both the quartic and quadratic divergences. 

To compute the trace (|o|) we need to be able to construct the three terms (a\e~ eH <> H\a) , (a\e~ eH ° Hq 1 H 2 \a) , 
and (a\e~ eH ° Ho\a) , where the states |a) are the elements of some basis for field perturbations. A moments thought 
reveals that the last two terms are simple to compute, since we have an explicit form for the operator H 2 , and the 
action Ho on any perturbation is trivial ((2^). The difficulty lies with the first term, the trace of H over its positive 
frequency subspace. This trouble arises because one does not have the H operator in a calculationally useful form. 



A. The H' 2 operator 



In order to calculate the mass correction, (fLO)) , we need the H 2 operator, which can be extracted from the classical 
perturbation dynamics. To second order in e a , the perturbation Lagrangian is 

C(e, de) = d t e a I ab d t e b - d z e a (V oWj -0j + W ab . ti ) e b - e a M ab e b , (11) 

where the tensors I(x),V(x),W(x), M(x) are functions of the static soliton fields. The classical time evolution of 
e a {x,t) is 

d 2 e a = [r 1 ] a " {Vbo-ijdij + <),\)„ : ,J), + diW bc . ti - M bc ) e c (12) 

= -H 2 oe. (13) 

And so the classical equations of motion for the perturbation yield an operator, H 2 , whose eigenvalues are the squares 
of the normal mode frequencies of perturbations about the soliton. 



B. Low normal modes of the soliton 



We have H 2 ; unfortunately its positive square root H cannot be easily represented. So, one cannot calculate all 
terms in the trace ( |l0| ) directly. However, let us note several points. First, the main contribution in (|l^) comes from 
the lowest frequency normal modes of the soliton; these are the modes where the difference between the (n\H\n) 
and (n\Ho\n) are greatest. Indeed several authors, e.g. ||, ||, go so far as to include only the contribution of 
the soliton zero modes to the trace. This is not an accurate approximation, however it is important to accurately 
compute the contributions of the lowest modes. Second, H 2 shares eigenvectors with H . So, although we cannot 
directly construct H, we can find its eigenvectors and eigenvalues if we can extract them from H 2 . Finally, it is clear 
that very short wavelength modes should be accurately described by the WKB approximation - their contribution to 
the renormalized mass correction should be small, local and analytically calculable via a derivative expansion. 

Our strategy is summarized in the following formula: 

Tr [O] = Tr„ [O] + Tr^ w [(1 - P)0] + Tt h>hm [(1 - P)0] , (14) 
where O = H^{H - H ) 2 . 

Here {|n)} is a small set of precise, low frequency normal modes for the soliton, computed numerically as explained 
below, and used to take the first trace. The remaining part of the trace is computed using a simple plane wave basis 
e jk x , corrected for overcounting using a projection operator 1 — P with P = ^2 n \n)(n\. The contribution from low 
momentum modes k < fc max is computed by straightforward matrix diagonalization. The final term from modes with 
k > fcmax is computed analytically via a local derivative expansion as explained below. 

The lowest frequency normal modes of H 2 may be extracted using the real-time method of |Q , or by using H 2 to 
drive a diffusion equation. The latter method has the advantage that for a limited set of modes the result converges 
more rapidly with integration time. For a fluctuation |e) about the static soliton, write 

d T \e)=-H 2 \e) . (15) 

Now, create an initial perturbation |e(0)). This may be written in terms of eigenvectors \n) of H 2 , with H 2 \n) = w 2 \n) . 

|e(0))=£e„|n) , (16) 



3 



so the solution of the diffusion equation (Jl5|) is 



|e(r))=^e-^e n ]n) . (17) 

n 

Starting in any state |e), and evolving it forward according to the diffusion equation (|l5|), the resulting state rapidly 
becomes dominated by the lowest frequency mode. 

To find the lowest frequency normal modes of H, one starts with an arbitrary perturbation field |e) and evolves 
it forward according to the diffusion equation, until the resulting field is an acceptably pure eigenvector of H 2 . Call 
this eigenvector |1) and store it. Next, take another perturbation, |e'), and project |1) out of it, and evolve it forward 
under the diffusion equation until one has another eigenvector of H 2 . Orthonormalize this with respect to |1), and 
call it |2), and store it. And so on: in this way one builds up an archive spanning the lowest frequency normal modes 
of perturbations about the soliton. Individual perturbations in this archive will be a mixture of nearby normal modes 
of H, since evolving the diffusion equation forward a finite length of time only separates out modes with substantially 
different frequencies. However, once created, this archive of wavefunctions will span the space of the the lowest 
frequency normal modes of H 2 , up to a cutoff frequency uj c < aj max , the highest frequency probed. From this mixed 
archive of low normal modes, one can create a near-perfect archive by diagonalizing the matrix (i\H creating a 
new archive of perturbations, out of linear combinations of the old ones. Any new perturbation field \i'), will be a 
near machine-perfect eigenvector of H 2 if its frequency satisfies 

where T is the finite time the diffusion equation evolved forward. All but the highest handful of modes extracted in 
this fashion will be near-perfect; those remaining will be discarded before further calculation. In the example shown 
later, we kept only the 9 lowest frequency modes (including the zero mode) as our set {\n}} of near-perfect low modes. 
The contribution of the near-perfect modes is easy to calculate, since 

H 2 \n)=u 2 n \n) =>■ H\n)=uj n \n) (19) 

and Hq acts on any perturbation via 

H \e) = J2 Vk 2 +m 2 \k)(k\e) . (20) 

k 

(Here the basis |fe) is an abbreviation for all possible plane wave perturbations of the field; for several component 
fields it carries an internal space index as well as a wavevector.) And so, in the basis of normal modes of H, it is 
straightforward to get each term in the trace. 



III. CALCULATING THE TRACE OVER ALL FLUCTUATIONS 



We wish to calculate the trace ( |Io| ) over all modes, however it is not practical to compute all of them as above. 
Instead we proceed by supplementing our near-perfect eigenmodes with a trace over plane waves. 
For the set of plane waves we simply construct the matrix 

H 2 kk , = {k'\H 2 \k) (21) 

for all plane waves k, k' satisfying k < fc max . We then diagonalize this matrix: 

H 2 kl =Oj d [H 2 D ] ddl O d , k , , (22) 

where O is an orthogonal matrix, and define the action of H via 

H kk ,=Ol d [Hl]l d ,O d , k , , (23) 
The eigenvectors of this matrix give an approximate basis of normal modes for the soliton: 

\d)= °^\ k ) ( 24 ) 



fc|<fc„ 
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and the eigenvalues [iff)] dd yield the squared frequencies of these approximate normal modes. These normal modes 
will be as perfect eigenvectors of H as can be constructed out of the limited region of Fourier space used. We may 
now straightforwardly compute each of the terms in the second trace in (|l4|). Most simply, if fc max is high enough, 
the trace will have converged and we are done. 

In order for this technique to provide an accurate mass correction, all terms in the trace, (|l0|), which contribute 
significantly to the mass correction must be included. That is, the subspace of normal modes spanned by this 
limited perturbation basis must cover all modes for which (n\(H — Hq) 2 / Ho\n) is large. It is possible to arrange 
this, because the finite contibutions from both far infrared and ultraviolet scattering modes are vanishing as their 
wavelengths become greatly larger or smaller than the soliton size. Let V^ i be the soliton volume scale, and let k so \ be 
a characteristic soliton wavevector scale, fc so i ~ V^^ d . infrared modes, with wavelength k < k m <C fc so i far away 
from the soliton give a total contribution ~ fc : d R {V so \/V\, ^). After removing analytically calculable corrections (see 
the next section), all ultraviolet modes with wavelength k > k vv 3> fc S oi give a total contribution ~ (l/fc uv ) 2 , if the 
Fourier transform of the soliton field dies away faster than 1/fc 2 for large k. (The solitons in which we are interested 
die exponentially with k.) So, the IR and UV cutoffs created by working in this truncated Fourier space create only 
controllable errors in the result. 



IV. EXAMPLE: QUANTUM CORRECTION TO THE MASS OF THE <^ 4 KINK 

As an example of this technique, we use it to calculate the first quantum mass correction to a well known example, 
the </> 4 kink in 1+1 dimensions. This model has been thoroughly explored; see pi for example. In the infinite volume 
limit the kink solution, and all its perturbation modes are known analytically |10(| . The calculation serves as a useful 
test since our numerical result can be compared with the exact analytic expression. Of course our technique can be 
straightforwardly be applied to any soliton bearing field theory in 1+1 dimensions. 

The Lagrangian for the 4 kink is 

C = \ {{d^f (d x 4>f) + \m 2 <? ~ ^ 4 ~ ^ , (25) 

with the usual double well potential in (f>, with minima at <f> = ±m/-\/A. The vacuum for this Lagrangian is (f> resting 
at one or the other minimum. The kink is a field 4>(x) which starts at one minimum at x — > — oo, and smoothly 
crosses over the potential hump at <f> = 0, to the other minimum as x — > oo. Minimizing the energy of such field 
configurations, one finds <fi s t{x) ( = m/vAtaiih(ma;/\/2), for the continuum theory, but the precise form of the 
function is unimportant to us.) Expanding in small perturbations about the vacuum, cf>(x,t) — <p vac + e(x,t), we 
obtain the Hq operator: 

H 2 = -dl + (V2m) 2 . (26) 
Similarly, expanding in small perturbations about the kink, we get 

tf 2 = -d 2 + (3A0 s 2 t (x)-m 2 ) . (27) 
which is the Schrodinger operator for a sech 2 potential. In order to calculate the mass correction, we discretise 



the system and and place it in a finite box. Equations (g6(),(g^) are discretised by the simple replacement d 2 e — > 
(e-t-i — 2ei + £i+i) /A 2 . Boundary conditions are handled by inserting a discontinuous jump of 2m/ V\ (the distance 
between the two vacua) at the box boundary. This choice of boundary deforms the soliton but the deformation 
vanishes exponentially with box size. The perturbations obey periodic boundary conditions. 

We begin by creating a stable soliton on the grid by relaxing kink-bearing initial conditions to a minimum energy 
configuration. The relaxed kink is shown in figure |l|. The field differs from its continuum value by no more than 
one part in 10 3 at any point in the box. Once we have the static soliton, we compute the lowest normal modes as 
described above, solving the diffusion equation by the Crank-Nicholson technique pi] . After diagonalization of the 
(i\H 2 \j) matrix, we obtained an archive of the lowest normal modes of perturbations about the soliton, each a perfect 
eigenvector of the discrete H 2 operator to machine accuracy. A selection of these normal modes is shown in figure |^. 

Using these near-exact low normal modes, together with the set of approximate normal modes constructed from 



the Fourier basis as in section III, we calculate terms in the mass correction trace, (JLOj) , using identities (|19|),(|20|), (|23j). 
Including all Fourier modes in the box, we obtained a total mass correction 8m = —0.9468. This is to be compared 
with the exact result (Rajaraman ||) of 8m — m/6 ^3/2 — 18/(n\/2)^J = —0.9422 for these values of parameters. 
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FIG. 1. The discretized 4> 4 kink. This, and all other figures, were calculated with m — 2, A = 1, on a 192 point grid, with a 
box size of 16 kink spatial units. 
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Thus this technique generates the right mass correction to about 0.5%, as good as could be hoped from a ~ 200 
point grid. As shown next, we do not need to include so many Fourier modes, since the high k contribution may be 
computed analytically. 



V. ANALYTIC COMPUTATION OF HIGH MOMENTUM CONTRIBUTION 



The contribution of high k Fourier modes to the mass correction can be calculated analytically using an expansion 
in k^ 1 . The result can then be used to improve the numerical result obtained using a limited range in Fourier space. 
The one dimensional theory studied here provides a test case for this approach. As we shall see, applying the lowest 
order correction makes the convergence with Fourier mode cutoff very rapid. 

As discussed above, the only nontrivial term is hTiH. We compute this using a heat kernel expansion as follows. 
First we set 

H 2 = -dl + U{x) + M 2 (28) 

with M = y/2m, the meson mass, and U — d 2 V/d<fi 2 — M 2 . The potential U defined this way vanishes at infinity. 
Next we write 

Performing the trace by summing over a complete set of plane wave states e lkx gives 

dk f , f dt , d, -t[-(ik+d*) 2 +u+M 2 ] 



where we used d x e tkx = d x + ik. We now separate the exponent into a 'free' part k 2 + M 2 and an 'interaction' 
part / = (—2ikd x — d 2 + U). We then expand e~ n out in powers of /, moving derivatives to the right noting that 
when a derivative reaches the extreme right, it gives zero. Thus we replace I — > U(x), I 2 — > U 2 — U" — 2ikd x U, 
I 3 — > (—Ak 2 )U" plus lower order terms in k. We then integrate over t, each power of t giving a factor of (k 2 + M 2 ) -1 . 
Keeping terms up to k~ 3 we find 



1 U 1 U 2 - U" 1 k 2 U h 



TrH = j-j d X Vk 2 + AP (I + - - &WTJ ^ - lWTJ ^ + -) (3D 

where we have ignored terms odd in k since they cancel if we adopt a symmetric cutoff k = \k\. 

Substituting ( |3l|) into the Cahill et al. formula (9), the first two terms of (^l|) are cancelled by the last two terms 
in (9). The remainder gives the high k contribution 

i r°° rik i r + °° ^/i/f 3 

thL m 1+u " ] =-^ k + ° ( > 1 >- (32) 

We may now use this formula to correct numerical results obtained with limited coverage in k space. That is, we 
compute the last term in (^) analytically as above. In this calculation we may ignore the term involving the projection 
operator P - as we shall argue, it vanishes exponentially as fc max is increased. To see this, note that the manipulations 
used in deriving ( |32j ) are actually still valid in the presence of the projector P, but with the introduction of a factor 
l(^l n )| 2 under the k integral. This involves sum of the squares of the Fourier transforms of the mode functions 
\n). But since the latter are smooth, the factor falls exponentially with fc max . 

Figure 3 shows the convergence of the 'improved' result as a function of the cutoff k. The horizontal line shows the 
full result obtained by including all modes in the box. The squares show the 'improved' result obtained for limited 
numbers of k vectors, and N labels the cutoff mode number, with k = 2tt N/L. The 'improved' result is within one 
per cent of the full result when we include modes up to N=7, or k ~ 5. 

Our technique clearly works very well for the (f> 4 kink in one spatial dimension, and could be applied with similar 
ease to any other one dimensional scalar field theory with kinks. The convergence of the 'improved' formula for 5m 
with the cutoff in k space augurs well for the prospect of performing analogous calculations in three dimensions. Just 
as here, the high k contributions can be computed analytically. And if the convergence is similar, we can expect 
to obtain accurate results for the soliton mass corrections using of order TV 3 ~ 10 3 modes in k space, which should 
certainly be feasible. 
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FIG. 3. Convergence of the computation of Sm with cutoff k = 2irN/L in a plane wave basis. Solid line shows the result 
including all modes in the box: dashed line indicates one per cent error, which is reached by N = 7. 



VI. CONCLUSION 



We have developed a technique for calculating the first quantum correction to the mass of a static soliton. This 
technique is implementable for arbitrary stable solitons within a broad class of Lagrangians. As an example, we 
applied this technique to the (j) 4 kink, where the result is known analytically: the numerical answer agreed well with 
the correct one. We hope to use this technique to calculate quantum corrections to the masses of more complicated 
solitons in three spatial dimensions. 
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